Precision Trial Designer showcase panel

European Institute of Oncology, Milan

European Institute of Oncology, Milan

In this analysis we will simulate a basket clinical trial run on patients affected by 3 diseases: breast, lung and gastric cancer. Samples will be collected and sequencing with Next-generation sequencing technology. Patients will be assigned to treatment groups based on the genomic characterization. The panel design is based on a recent publication that defines the guidelines for precision medicine and the genes that should be included into diagnostic NGS assays for the previously indicated tumor type.

Consensus on precision medicine for metastatic cancers: a report from the MAP conference C.

BASKET TRIAL PANEL SHOWCASE

OBJECTIVES

The objective of this analysis is to:

1. Simulate the panel “detection power” (defined as the number of patients with at least one alteration detectable by the genomic panel) in different tumor types.
2. Simulate a genomic driven clinical trial: As we might risk to have “unbalanced” group treatment, we would like to know in advance the patient allocation frequency to different treatment arms based on the alteration frequency in the cancer population;
3. Estimate the number of patients needed to be screened;
4. Optimize the size of the panel and therefore reduce experimental cost;
5. As last, illustrate the usage of PTD in a concrete showcase context.

Initial assumption

The panel is designed based on the guidelines indicated in the recent publication:

Consensus on precision medicine for metastatic cancers: a report from the MAP conference C.

We had to make some adjustments and interpretations in order to fit the guidelines into a panel design that would be input by the PTD package.

Gene name adoption

PTD requires gene names to be written in HGNC format. Since the software will download data from different web databases, we need to make sure that gene names follow the international convention. While scanning the paper, we noticed that the following genes were encoded in the paper with names that do not match the HGNC standards. We have therefore apply the following interpretations and corrections:

Paper annotation	Conversion for the analysis
NTKR (Possible typo?)	NTRK1, NTRK2, NTRK3
NTRK	NTRK1, NTRK2, NTRK3

Moreover, the paper states that for gastric tumor, MSI should be tested. We can infer this information testing the following genes: MLH3, MLH1, MSX1, MSX2

Gene alterations

In this analysis we will only consider the following alteration:

• SNV
• CNA
• Fusions (only for lung tumor since they were selected to be lung cancer specific alteration)

Changes in expression are not included in this analysis.

Drug Treatment

For demonstrative purposes we are simulating a clinical trial by associating altered genes to “general drug families” as a possible treatment.

Tumor types

The paper we are basing our work on, defines the guidelines for 3 diseases:

• Breast cancer;
• Lung cancer;
• Gastric Cancer;

We need to convert these diseases into tumor ID compatible with the TCGA standard. The conversion is applied as follows:

Publication definition	TCGA tumor ID conversion adopted
Breast Cancer	brca (Breast invasive carcinoma)
Lung Cancer	luad (Lung adenocarcinoma)
Lung Cancer	lusc (Lung squamous cell carcinoma)
Gastric Cancer	stad (Stomach adenocarcinoma)
Gastric Cancer	coadread (colon)

We are making the assumption that the cancer population of reference is composed entirely and uniquely by this cancer subtypes

# Define tumor types 
tumor_types <- c("brca", "luad", "lusc", "stad", "coadread")

Sampling Populations

In this section, we will attempt to estimate the sample size of each cancer type in a population of reference. For this simulation we will use a population of reference sampled according to the statistics published by seer.cancer.gov

• http://seer.cancer.gov/csr/1975_2013/results_single/sect_a_table.01_2pgs.pdf
• local copy of the file extra_docs/sect_a_table.01_2pgs.pdf

Seer reports that, between 2009 and 2013 the tumor population amounted to:

• Breast cancer: 306665 people;
• Lung cancer: 256350 people;
• Gastric cancer:
  • Stomach cancer: 33174 people;
  • Colon and rectum cancer: 185758 people;

Breast carcinoma

Frequency of the population is adjusted by Breast carcinoma according to the following table sect_04_table.25_2pgs.pdf:

• brca: 99.4%

Lung subtypes

Frequency of the population is adjusted by lung tumor subtypes according to online link, local file:

• Lung Adenocarcinoma: 45% of lung population
• Lung squamous and transitional cells: 23.0% of lung population

Gastric subtypes

Frequency of the population is adjusted by stomach tumor subtypes according to online link, local copy

• Stad: 82%

Frequency of the colon cancer are adjusted according to online link, local file

• Coadread: 92.8%

We will therefore use these statistics as our population of reference using the following approximations:

set.seed(7)

# assign number of patients
luad_tot = ceiling(256350*0.45) #adjustement coefficient specific to the % of subtype tumor in a lung cancer population
lusc_tot = ceiling(256350*0.23) #adjustement coefficient 
stad_tot = ceiling(33174*0.82) #adjustement coefficient 
coadread_tot = ceiling(185758*0.928) #adjustement coefficient 
brca_tot = ceiling(306665 * 0.994)

# total
all_tot = brca_tot+luad_tot+lusc_tot+stad_tot+coadread_tot

# calculate frequencies
brca_freq = brca_tot/all_tot
luad_freq = luad_tot/all_tot
lusc_freq = lusc_tot/all_tot
stad_freq = stad_tot/all_tot
coadread_freq = coadread_tot/all_tot

TCGA tumor ID	population considered	frequency
brca	304826	0.449111
luad	115358	0.169961
lusc	58961	0.0868693
stad	27203	0.0400791
coadread	172384	0.2539795

ANALYSIS

In the following phases we will proceed by designing the panel and visualizing its performance through different types of plots.

Design phase

First step is to create the panel design using an excel file according to the standards of the package. The panel design is loaded into R and can be previewed in the table below.

library(DT)
library(knitr)
library(ggplot2)
library(PrecisionTrialDesigner)

# Helper functions
`%++%` <- function(a,b) paste0(c(a,b), collapse = " ")

# calculate the right confidence interval of a proportion
rightCIprop <- function(p, sample_size){
  se <- sqrt( p*(1 - p)/sample_size)
  right <- p + 1.96*se + 0.5/sample_size
  if (right > 1) right <- 1
  return(right)
}

# calculate the left confidence interval of a proportion
leftCIprop <- function(p, sample_size){
  se <- sqrt( p*(1 - p)/sample_size)
  left <- p - 1.96*se - 0.5/sample_size
  if (left < 0) left <- 0
  return(left)
}

# Import excel file in R
panel_design <- gdata::read.xls("../Tables/tb1_PTD_panel_showcase.xlsx")
# preview the input table
datatable(panel_design)

We can now run the main function of the R package and:

1. Import the panel design into a CancerPanel object;
2. download Alteration information for the drug types of interest;
3. Adapt the alteration information to the CancerPanel design;

# 1. Create the panel Object using PTD package - DESIGN PHASE
PTD_panel <- newCancerPanel(panel_design)
# 2. DOWNLOAD PHASE
PTD_panel <- getAlterations(PTD_panel, tumor_type=tumor_types)
# 3. SIMULATION PHASE
PTD_panel <- subsetAlterations(PTD_panel)

#save a backup copy to speed up future analysis
saveRDS(PTD_panel, file="../Temp/PTD_panel.rds")

At this point PTD_panel contains all the information we need to visualize the panel performance.

Filter non pathogenic mutations

Remove NON pathogenic mutation fetched from cosmic dataset (01-05-2017).

Cosmic mutations where downloaded from the cosmic dataset on the 18-05-2017.

# save the number of entries in the TCGA beofere any filtering
TCGA_original_number <- nrow(unique(PTD_panel@dataFull$mutations$data))

# Clean TCGA popualation running a cosmic filter
cosmic <- readr::read_csv("../external_resources/cosmic_df_PTD.csv")

## Warning: Missing column names filled in: 'X1' [1]

## Parsed with column specification:
## cols(
##   X1 = col_integer(),
##   Position = col_integer(),
##   CDS.Mutation = col_character(),
##   AA.Mutation = col_character(),
##   `Mutation ID (COSM)` = col_integer(),
##   Count = col_integer(),
##   Type = col_character(),
##   gene = col_character()
## )

cosmic_df <- data.frame(gene_symbol = cosmic$gene
                        , amino_acid_change = sub("^p." , "" , cosmic$AA.Mutation)
                        , stringsAsFactors = FALSE)
PTD_panel <- filterMutations(PTD_panel , filtered = cosmic_df , mode = "keep")

## Subsetting fusions...

## Subsetting mutations...

## Subsetting copynumber...

• Filtering out ** 33.5 % ** of the TCGA original SNVs.

Exclude predefined genomic regions from the simulation

A useful feature of PTD is the possibility to remove with precision certain genomic regions from the study simulation. This is particularly relevant in the clinical routing as, we might be interested for example only in a subset of mutation in a particular portion of a gene. For example, EGFR T790M and PIK3CA exon 20 mutations are associated with resistance to 1st generation EGFR inhibitors. In this analysis we will model the SNV mutations in PIK3CA excluding SNV located in exon 20, as well as exclusding EGFR T790M.

This can be achieved through the shiny app web interfance using the function panelOptimizer or programmatically, by applying a genomic filter:

# PIK3CA, Exon20 on the hg19 genome assembly is mapped at the 
# coordinates: chr3: 178,866,311-178,957,881
# ==============================================================================
# Let's create a bed file with the genomic region to exclude
bed <- data.frame(chr = "chr3"
                  , start = 178948013
                  , end = 178948164
                  , stringsAsFactors = FALSE)
knitr::kable(bed , row.names = FALSE)

chr	start	end
chr3	178948013	178948164

PTD_panel <- filterMutations(PTD_panel , bed = bed , mode = "exclude")

## Subsetting fusions...

## Subsetting mutations...

## Subsetting copynumber...

# EGFR T790M can be removed by specifing the mutations using the a.a. format
# ==============================================================================
toberemoved <- data.frame(gene_symbol = "EGFR"
                          , amino_acid_change = "T790M"
                            )
PTD_panel <- filterMutations(PTD_panel , filtered =  toberemoved, mode = "exclude")

## Subsetting fusions...

## Subsetting mutations...

## Subsetting copynumber...

# for help and more details see ?filterMutations

Adjust CNA threshold to be clinically relevant

CNA data provided by cBioPortal reports the final step in CNV calls from GISTIC pipeline (CIT). CNVs are reported as either

• -2: deletion
• -1: shallow deletion
• 0: normal
• 1: shallow amplification
• 2: amplification

Values of -2 or 2 are considered alterations, while the rest is considered as normal copies. For research purposes, these criteria are very useful to discover new potential gene copies alterations but a clinical setting requires stronger evidence for a call, especially in case of amplification. The SHIVA trial protocol, for instance, in order to consider clinically relevant a AKT amplification, it requires at least the presence of 5-6 copies. A similar granularity is not possible with the final GISTIC results but it can be obtained up to a certain degree using the original data provided by the GDAC Firehose (the same used by cBioPortal).

We downloaded GISTIC run 2016_01_28 all_thresholded.by_genes.txt files for 32 tumor types and normalized the segmented value for each patient by their respective thresholds. Each sample of the TCGA cohort has a deletion and amplification call threshold calculated according to estimated ploidy and purity of the original sample. For example, a low purity sample will have lower threshold for calling a CNV compared to a 100% purity sample because of the normal diploid cells contamination. These thresholds represents the points were an amplification (2) or deletion (-2) are called.

This normalization reformats the data as absolute number of copies in the following way:

• 0-1: deletion
• 1.1: shallow deletion
• 2: normal
• 2.9: shallow amplification
• 3-Inf: amplification

While data from 1.1 to 2.9 are the same as before (considered normal on default), now 1 is the point of -2 deletion and 3 the point of +2 amplification. We have now the full spectrum of amplification from > 3 that can be used to fine-tune the clinically relevant threshold (4, 5 ,6 copies for example).

# SETTINGS
# ==============================================================================
# DEFINE GENES: define genes we are interested for the CNA analysis
CNA_genes <- PTD_panel@arguments$panel %>%
  dplyr::filter(alteration == "CNA") %>%
  dplyr::select(gene_symbol) %>%
  .[[1]]

# PATH: set global variable with the path to.rds files with CNA raw data
PATH_2CNA_FILES <- "/Users/ale/IEO/projects/IEO_research/R_package_panel_and_projection/functionfactory/cna_firehose_rds/"
#PATH_2CNA_FILES <- "~/Dropbox/git_guida/PrescriptionTrialDesigner_publication/version_6.0/Temp/cna_firehose_rds"

# helper function that reads in a firehose RDS file, filters by gene, and apply the threshold configured
# ==============================================================================
get_cna_long <- function(tumor_types, PATH_2CNA_FILES, genes_chr, del_thr, amp_thr){
  
  
  # CHECKS
  # --------------------------------------------------------------------------------
  # check that PATH actually leads to a dir
  if(!dir.exists(PATH_2CNA_FILES)){ stop("path to a dir that cannot be found")}

  # Check that all the files are found for each tumor type required in the analysis
  input_lgl <- purrr::map_lgl(tumor_types
                              , ~ file.exists(file.path(PATH_2CNA_FILES
                                                     , paste0(toupper(.) , "_PTD.rds")
                                                     )
                                              )
                              )
  
  if(!all(input_lgl)){
    stop(paste("Files for some of the tumor types where not found:"
               , tumor_types[-which(input_lgl)])
         )
  }
 
  # Checks
  if(!is.character(tumor_types) | is.null(tumor_types) | length(tumor_types) ==0){
    stop(paste("Tumor type is not confirm", tumor_types))
  }


  # Define an internal function to do the job
  # -------------------------------------------------------------------------------
  get_CNA_long_onetumor <- function(one_tumor, PATH_2CNA_FILES, genes_chr, del_thr, amp_thr){
    
    # read in CNA file
    # myfile <- readRDS(paste0(PATH_2CNA_FILES, toupper(one_tumor) , "_PTD.rds"))
    # if(!file.exists(myfile)){
    #   message(myfile)
    #   return(NULL)
    # }
    # mat <- tryCatch(readRDS(myfile) , error = function(e) stop(myfile))
    
    mat <- readRDS(file.path(PATH_2CNA_FILES, paste0(toupper(one_tumor) , "_PTD.rds")))
    
    # extract Patients
    pats <- colnames(mat)

    # Filter for gene of interest
    mat <- mat[ genes_chr , ]

    # Melt and adjust
    output <- reshape2::melt(mat , value.name = "CNAvalue")
    colnames(output) <- c("gene_symbol" , "case_id" , "CNAvalue")
      output$tumor_type <- one_tumor
      
    # Choose cutoff values
    # -----------------------------------------------------------------------------
    # With default threshold data are now like this: 
        # [0 , 1.1) "deletion"
        # 1.1 "shallow deletion"
        # 2 "normal"
        # 2.9 "shallow amplification"
        # (2.9 , Inf] "amplification"
    output$CNA <- ifelse( output$CNAvalue <del_thr 
                    , "deletion" 
                    , ifelse(output$CNAvalue > amp_thr , "amplification" , "normal"))
    # select
    output <- output %>%
      dplyr::select(gene_symbol , CNA , case_id , tumor_type, CNAvalue) %>%
      dplyr::mutate(gene_symbol = as.character(gene_symbol)) %>%
      dplyr::mutate(case_id = as.character(case_id)) 
    
    return(output)
  }

  # Main body of the function
  # -----------------------------------------------------------------------------
  cna_long_df <- purrr::map(tumor_types
                             , ~ get_CNA_long_onetumor(.
                                                       , PATH_2CNA_FILES
                                                       , genes_chr
                                                       , del_thr
                                                       , amp_thr
                                                       )
                             ) %>%
    dplyr::bind_rows()
    
  return(cna_long_df)  

}
  

# RUN MAIN FUNCTION
# ==============================================================================
cna_long_df <- get_cna_long(tumor_types
                            , PATH_2CNA_FILES
                            , CNA_genes
                            , del_thr = 1.1
                            , amp_thr = 4
                            )


# Create Sample
# ----------------------
# get list of tumors
tumor_chr <-  cna_long_df$tumor_type %>% 
  unique %>% 
  as.character

# prepare funtion to extract sample names for each tumor
sampleXtumor <- function(x, df){
  df %>% 
    dplyr::select(case_id, tumor_type) %>%
    dplyr::filter(tumor_type %in% x) %>% .[["case_id"]] %>%
    as.character() %>%
    unique()
}
# Extract samples for each tumors
Samples <- purrr::map(tumor_chr, ~ sampleXtumor(., df=cna_long_df))
names(Samples) <- tumor_chr

# Check if there are problems. 
testthat::expect_equal(length(unlist(Samples)), length(unique(cna_long_df$case_id)))

# Create list element for CNA for the DataFull argument of the panel (panel@DataFull)
# ==============================================================================
copynumber <- list(  data = data.frame(cna_long_df)
                   , Samples = Samples)


# Preview changes before committing
PTD_cna_pop_comparison <- purrr::map2( 
    .x = list( PTD_panel@dataFull$copynumber$Samples, copynumber$Samples )
  , .y = list( "TCGA", "Firehose" )
  , .f = ~ purrr::map_int(.x, length) %>%
      data.frame(n_patients = .) %>%
      tibble::rownames_to_column("tumor_type") %>%
      dplyr::mutate(source = .y)
  ) %>% 
  dplyr::bind_rows() %>%
  ggplot(aes(tumor_type, n_patients, fill=source)) + 
    geom_bar(stat="identity", position= position_dodge()) +
    coord_flip() + 
    ggtitle("number of patients in the CNA datasets")

# Preview
PTD_cna_pop_comparison

apply changes to the TCGA population object.

# Apply changes
PTD_panel@dataFull$copynumber <- copynumber
PTD_panel <- subsetAlterations(PTD_panel)

## Subsetting fusions...

## Subsetting mutations...

## Subsetting copynumber...

Visualisation and Simulations of the panel performance

In this section we will start running some of the visualization functions to investigate how well the panel would perform in different conditions. Before starting out investigation we want to visualize the sample number by tumor type to keep in consideration any sampling bias due to the population of origin. This bias can be corrected using the option tumor.wieights or tumor.freqs.

The following table shows the number of samples (y-axis) collected for this analysis for each tumor type considered.

coverageStackPlot(PTD_panel
                  , alterationType=c("mutations" , "copynumber") 
                  , var="tumor_type", noPlot=TRUE)$Samples %>%
  data.frame(Samples=.) %>%
  tibble::rownames_to_column("tumor_type") %>%
  knitr::kable()

tumor_type	Samples
brca	976
coadread	220
luad	230
lusc	178
stad	393
all_tumors	1997

Let’s export the data and make a summary table to resume differences in populations:

# get data for table
tcga_pop <- coverageStackPlot(PTD_panel
                  , alterationType=c("mutations" , "copynumber") 
                  , var="tumor_type"
                  , noPlot=TRUE)

# calculate frequencies
tcga_freqs <- tcga_pop$Samples[-length(tcga_pop$Samples)]  / tcga_pop$Samples["all_tumors"]

# make it into a df
population_df <- data.frame(tumor=c("brca", "luad", "lusc", "stad", "coadread")
          , seer=c(brca_freq, luad_freq, lusc_freq, stad_freq, coadread_freq)
          , tcga = c(tcga_freqs["brca"], tcga_freqs["luad"], tcga_freqs["lusc"], tcga_freqs["stad"], tcga_freqs["coadread"])
)

# barplot, with only the genes in the panel
population_df %>%
  tidyr::gather(Population, Frequency, seer:tcga) %>%
  ggplot(aes(x=Population, y=Frequency, fill=tumor)) + 
  geom_bar(stat="identity") + 
  facet_grid(tumor ~ ., space="free", scales="free") +
  ggtitle("Comparison of the initial population composition") + 
  coord_flip()

library(svglite)
ggsave(filename="../Figures/fig2A.svg", plot=last_plot(), device = "svg")

## Saving 10 x 8 in image

Comparison table between SEER initial population and TCGA initial population

TCGA tumor ID	Seer population considered	Seer frequency 1	TCGA initial population	TCGA sample frequency 2
brca	3049	0.449111	976	0.4887331
luad	1154	0.169961	230	0.1151728
lusc	590	0.0868693	178	0.0891337
stad	273	0.0400791	393	0.1967952
coadread	1724	0.2539795	220	0.1101652

Fusion detection in Lung cancer

The previous plot does not include fusions, as the fusion in the panel are specific to the lung tumor. If we want to visualize the number of samples with available information for mutations, copynumber and fusions (each patients will have information for each category).

coverageStackPlot(PTD_panel
                  , alterationType=c("mutations" , "copynumber", "fusions") 
                  , var="tumor_type"
                  , tumor_type = c("luad", "lusc")
                  )

Panel detection power in the overall population

One common question that is often considered when trying to design custom target enrichment panel, is the percentage of the population that would be “captured” by the design of the panel. The following plot shows the number of alteration (x-axis) found in the sample (patient) population (y-axis). The panel detection power can be defined as the number of patients with at least one alteration detectable by the genomic panel. In this simulation we are interested in visualizing the detection power in a population background composed by patients from breast, lung, and gastric cancer. For simplicity fusion information are not included in the analysis.

coveragePlot(PTD_panel
             , alterationType=c("mutations" , "copynumber")
             # adjust population distribution by applying a weighted sampling bootstrap
             , tumor.freqs = c(brca=brca_freq
                               , luad=luad_freq
                               , lusc=lusc_freq
                               , stad = stad_freq
                               , coadread = coadread_freq
                               )
            )

Simulate the panel performance with a weighted bootstrap simulation run 200 times.

The previous plot shows the percentage of patients (y-axis) found to have at least one, two, three, etc, alterations (x-axis). The current design is capable of detecting 0.65% of the population of reference.

Panel detection-power in a tumor type specific population background

In the following plot we want to identify the performance of our panel for each single tumor types starting from out population of reference composed by the 3 cancers groups.

coveragePlot(PTD_panel
             , alterationType=c("mutations" , "copynumber")
             , grouping = "tumor_type"
             , tumor.weights = c(brca=brca_tot
                               , luad=luad_tot
                               , lusc=lusc_tot
                               , stad = stad_tot
                               , coadread = coadread_tot
                               )
             )

To allow and easy comparison of the detection power in different tumor types we have redesign the plots into one unique barplot.

# get data from the panel simulation
dt_power_lt <- coverageStackPlot(PTD_panel
                  , alterationType=c("mutations" , "copynumber") 
                  , var="tumor_type"
                  , noPlot = TRUE
              )

# convert to dataframe
dt_power_lt$plottedTable <- data.frame(dt_power_lt$plottedTable)

# calculate frequency
df <- data.frame()
for(i in 1:length(dt_power_lt$plottedTable)){
  
  # find tot. number of samples in matching tumor
  index <- which(names(dt_power_lt$Samples) %in% names(dt_power_lt$plottedTable[i]))
  
  if(length(index) != 0 ){
    # calculate frequency
    # freq <- dt_power_lt$plottedTable[i] / dt_power_lt$Samples[index]
    freq <- dt_power_lt$plottedTable[i] / dt_power_lt$Samples["all_tumors"]
  }
  # export
  df <- rbind(df
              , data.frame( tumor_id=names(dt_power_lt$Samples[index])
                          , freq=as.numeric(freq)
                          , n_samples = dt_power_lt$Samples[index]
                           )
              )
}

# calculate confidence intervals
df$leftCI <- purrr::map2_dbl(df$freq, df$n_samples, .f=leftCIprop)
df$rightCI <- purrr::map2_dbl(df$freq, df$n_samples, .f=rightCIprop)

# sort dataframe
df <- df %>% dplyr::arrange(freq) 

# plot
df %>%
  ggplot(aes(x=reorder(tumor_id, freq), y=freq*100)) + 
  geom_bar(stat="identity") + 
  scale_y_continuous(limits = c(0,50)) +
  geom_errorbar(aes(   ymax= leftCI*100
                     , ymin= rightCI*100)
                 , position=position_dodge(0.9)
                 , size=0.2
  ) +
  geom_text(aes(x=df$tumor_id
                , y=df$freq*100+0.15
                , label = sapply(df$freq*100
                             , function(x) paste(c(round(x, digits=2), "%")
                                                 , collapse=" ")
                            )
                  ) , position = position_dodge(width = 1)
                    , vjust = 0.5
                    , hjust = -1.3
              ) + 
  coord_flip() + 
  labs(title="Genomic panel detection-power in different tumor types.") +
  ylab(label = "detection power indicated as the overall population percentage") +
  xlab(label = "tumor type (TCGA ID)")

Patient treatment allocation

We can also easily address the issue of the allocation frequency in the different treatment groups. The following plot shows the number of patients with at least one, two or more alterations that can be targeted by each drug in the trial. This is particularly important in clinical trials when you want to predict the maximum number of patients that could be allocated into a drug treatment and identify potential cohorts with few patients allocated.

coveragePlot(PTD_panel
             , alterationType=c("mutations" , "copynumber")
             , grouping="drug"
             # adjust population distribution by applying a weighted sampling bootstrap
             , tumor.freqs = c(brca=brca_freq
                               , luad=luad_freq
                               , lusc=lusc_freq
                               , stad = stad_freq
                               , coadread = coadread_freq
                               )
             , colNum = 2 # number of columns in the plot
             , cex.main = 1 # size of the plot title
)

To make it simpler to compare, we will plot the detection power for each drug type in a unique plot. The following code will put all the previous plot together in one unique barplot.

# get data from the panel simulation
drug_df <- coveragePlot(PTD_panel
             , alterationType=c("mutations" , "copynumber")
             , grouping="drug"
             # adjust population distribution by applying a weighted sampling bootstrap
             , tumor.freqs = c(brca=brca_freq
                               , luad=luad_freq
                               , lusc=lusc_freq
                               , stad = stad_freq
                               , coadread = coadread_freq
                               )
             , colNum = 2 # number of columns in the plot
             , cex.main = 1 # size of the plot title
             , noPlot=TRUE
 )$plottedTable %>% data.frame()

# name first column with frequency of patients with at least one alteration that makes the patients elibigle for treatment with a specific drug
names(drug_df)[1] <- "freq"

# generate plot 
drug_df %>%
  tibble::rownames_to_column("drug_name") %>%
  dplyr::filter(!grepl("no_drug", drug_name)) %>%
  ggplot(aes(x=reorder(drug_name, freq), y=freq*100)) + 
  geom_bar(stat="identity") + 
  coord_flip() + 
  labs(title="Max frequencies of patients elibigle for each drug treatment") + 
  ylab(label = "detection power indicated as the overall population percentage") +
  xlab(label = "Drug name")

Generate the same plot this time running a weighted bootstrap simulation.

# Run simulation
# ------------------------------------------------------------------------------
N_SIMULATIONS = 200

simul <- lapply(1:N_SIMULATIONS , function(x) {
  drug_df <- coveragePlot(PTD_panel
             , alterationType=c("mutations" , "copynumber")
             , grouping="drug"
             # adjust population distribution by applying a weighted sampling bootstrap
                          , tumor.weights = c(brca=ceiling(brca_tot/100)
                                 , luad=ceiling(luad_tot/100)
                                 , lusc=ceiling(lusc_tot/100)
                                 , stad = ceiling(stad_tot/100)
                                 , coadread = ceiling(coadread_tot/100)
                                 )
             , colNum = 2 # number of columns in the plot
             , cex.main = 1 # size of the plot title
             , noPlot=TRUE
             , maxNumAlt = 1
 )
  
  return(drug_df$plottedTable/drug_df$Samples[1])
})
  
# Convert to dataframe
drug_vector <- unique(PTD_panel@arguments$panel$drug)

# Create a custom function to manage the problem with the function output CoveragaPlot not 
# outputting the same number of ites in the simualation  The function takes the labels, the 
# vector and rematches, adding NA when a match is missing
value2slot <- function(value_vector, labels_vector){

  # Cycle through the vector
  temp_vector <- rep(0, times=length(labels_vector))
  for (i in 1:length(labels_vector)){
    # for each position in the vector fetch
    value <- value_vector[i]
    label <- rownames(value_vector)[i] %>% gsub("drug: ", "", .)
      
    #get location index
    index <- labels_vector %in%  label
    # qc
    if(sum(index)>1) {stop("label mateched to more than one location")}
    # replace value to that index position
    temp_vector[which(index)] <- value
  }
  
  # add labels to the vector
  names(temp_vector) <- labels_vector
  
  return(temp_vector)
}
# test the function
#value2slot(simul[200][[1]], drug_vector)

# apply transformation to make it into a dataframe
simul_df <- purrr::map(simul, ~ value2slot(., drug_vector)) %>%
  data.frame %>%
  t()
# remove undesired rownames
row.names(simul_df) <- NULL

# GET STATS
# ------------------------------------------------------------------------------
# Calculate for each n° of alteration, the mean
detectionPower <- colMeans(simul_df)

# Calculate lower CI for bootstrap
leftCI <- apply(simul_df,2,function(x){left <- quantile(x, 0.025)})

# Calculate upper CI for bootstrap
rightCI <- apply(simul_df,2,function(x){right <- quantile(x, 0.975)})

# Preapre for column with # alterations
n_alterations <- 1:length(detectionPower)

# Put them all together in a data frame
dtpow <- data.frame(n_alterations, detectionPower, leftCI, rightCI) %>%
  tibble::rownames_to_column("drug_name") %>% 
  dplyr::filter(!drug_name %in% "no_drug") # remove no_drug from the dataset

# Generate the plot
PTD_detection_image <- dtpow %>%
  ggplot(aes(x=reorder(drug_name, detectionPower), y = detectionPower)) + 
  geom_bar(stat="identity" , position = position_dodge()) + 
  scale_y_continuous(limits = c(0,.5)) +
  #scale_x_discrete(labels = dtpow[order(dtpow$detectionPower, decreasing = TRUE),"drug_name"]) + 
  #scale_x_continuous(breaks = n_alterations, trans = "reverse") +
  geom_errorbar(aes(ymax=rightCI
                    , ymin=leftCI)
                , position=position_dodge(0.9)
                , size=0.2
  ) + 
  ggtitle("Max frequencies of patients elibigle for each drug treatment") +
  ylab(label = "Detection power indicated as the overall population percentage") +
  xlab(label = "Drug name") +
  annotate("text"
           , x=reorder(dtpow$drug_name, dtpow$detectionPower)
           , y=dtpow$detectionPower + 0.15
           , label= sapply(dtpow$detectionPower*100
                           , function(x) paste(c(round(x, digits = 2), "%")
                                               , collapse = " ")
           )
  ) +
  coord_flip()

# Show image
PTD_detection_image

Panel performance evaluation

The design phase of a panel requires the user to make the important decision about what genes to include. PTD provides a suite of tools envisioned to aid the user during the decision process and evaluate the performance of the panel.

Rank genes in the design in terms of genomic space/detection power

The following plot was studied to investigate the relationship between the “cost” of adding new genes to the panel (as a function of genomic space to be included among the genomic regions to be targeted) and the detection power in terms of the mean number of alterations found in the samples (cancer population), for that specific gene in a specific set of tumor types.

In order to optimize a panel design, we might be interested to know:

• What genes are most altered in the tumor types of interest;
• What genes are “cheaper” to be sequenced (in terms of genomic spaces, or also size of the targeted region);
• What genes have the highest number of (alteration frequency) / (genomic space).

All tumors together

The following saturation plot shows the overall performance of the panel in a population that includes all tumor types of interest (balanced by tumor frequency in the population of reference). Genes can be ranked based on the ratio (alteration frequency) / (genomic space). Genes with highest mean number of alteration per sample and smallest size are listed first. The following genes in the plot are added one by one incrementally. This means that on the y-axis, every time we add a new gene, we are adding only the incremental contribution of alteration for that specific gene. This was done to prevent from counting twice (or more than once), patients with genes that were already found altered.

saturationPlot(PTD_panel
            , alterationType=c("mutations" , "copynumber")
            , adding="gene_symbol"
            # adjust population distribution by applying a weighted sampling bootstrap
            , tumor.weights=c(brca=brca_tot
                             , luad = luad_tot
                             , lusc = lusc_tot
                             , stad = stad_tot
                             , coadread = coadread_tot
                             )
            , adding.order = "rate"
            , legend = "out"
            )

Each single cancer type will have different gene alteration frequencies, and therefore the order of the genes can change accordingly. In the following set of plots we will change the gene adding order, using asbolute ranking. This parameters, will not add genes based on their optimal ratio, but will priorities genes with the absolute highest number of alterations found in each single tumor group of interest

Breast Cancer

saturationPlot(PTD_panel
            , alterationType=c("mutations" , "copynumber")
            , adding="gene_symbol"
            , tumor_type = "brca"
            , adding.order = "absolute"
            )

Lung Cancer, mutations and copy number alterations

saturationPlot(PTD_panel
            , alterationType=c("mutations" , "copynumber")
            , adding="gene_symbol"
            , tumor_type = c("luad", "lusc")
            , adding.order = "absolute"
            )

Lung Cancer, only fusion alterations

saturationPlot(PTD_panel
            , alterationType=c("fusions")
            , adding="gene_symbol"
            , tumor_type = c("luad", "lusc")
            , adding.order = "absolute"
            )

## Warning in dataExtractor(object = object, alterationType =
## alterationType, : The following tumor types have no alteration to display:
## lusc

Gastric Cancer

saturationPlot(PTD_panel
            , alterationType=c("mutations" , "copynumber")
            , adding="gene_symbol"
            , tumor_type = c("stad", "coadread")
            , adding.order = "absolute"
            )

Rank drugs in the design in terms of genomic space/patient allocation

Similarly to the previous set of plots, we can investigate the cost of adding a new drug in the study in terms of extra genomic space to be targeted by the panel with respect to the alteration frequency of population that can be captured.

saturationPlot(PTD_panel
            , alterationType=c("mutations", "copynumber")
            , adding="drug"
            , y_measure="absolute"
            , tumor.weights = c(brca=brca_tot
                               , luad=luad_tot
                               , lusc=lusc_tot
                               , stad = stad_tot
                               , coadread = coadread_tot
                               )
             
            )

Similarly, we can separate the saturation line by tumor_types. (Note: In the following plot we do not need to apply the tumor.wight option, as we are grouping by tumor_types, Note2: no_drug indicates frequency associated with driver genes which are therefore not associated to any drug treatment in this specific design).

saturationPlot(PTD_panel
            , alterationType=c("mutations", "copynumber")
            , grouping="tumor_type"
            , adding="drug"
            , y_measure="absolute"
            , legend = "out"
            )

To improve the readability of the plot, we can look at each single tumor type individually.

Breast Cancer

# Generate plots for each tumor types
saturationPlot(PTD_panel
            , alterationType=c("mutations", "copynumber")
            , grouping="tumor_type"
            , adding="drug"
            , y_measure="absolute"
            , tumor_type = "brca"
) 

Lung Cancer

saturationPlot(PTD_panel
            , alterationType=c("mutations", "copynumber")
            , grouping="tumor_type"
            , adding="drug"
            , y_measure="absolute"
            , tumor_type = c("luad", "lusc")
)

Gastric Cancer

saturationPlot(PTD_panel
            , alterationType=c("mutations", "copynumber")
            , grouping="tumor_type"
            , adding="drug"
            , y_measure="absolute"
            , tumor_type = c("coadread", "stad")
) 

Panel detection power of each alteration type in the overall population

SNV and indels detection are today part of standard analysis for NGS application.However, CNA (Copy number alteration) detection, or fusion analysis are still not standardized analysis and they do not represent an easy implementation. In the following plot we will investigate the contribution of the detection power of each alteration type in case. This will help us quantify, in case we decide to opt out the detection of some of the alteration type, the influence of each single alteration type in the detection power.

saturationPlot(PTD_panel
            , alterationType=c("mutations" , "copynumber", "fusions")
            , adding="gene_symbol"
            , tumor.weights=c(brca=brca_tot
                             , luad = luad_tot
                             , lusc = lusc_tot
                             , stad = stad_tot
                             , coadread = coadread_tot
                             )
            ,grouping="alteration_id"
            )

## Warning in dataExtractor(object = object, alterationType =
## alterationType, : The following tumor types have no alteration to display:
## stad, coadread

Drug targeting frequency plot

In case some of the drugs are multitarget (in this panel design, we are using drug families), it might be interesting to investigate what are the most targeted genes in our population of reference. We can visualize this information as an independent stack plot for each tumor type.

par(mfrow=c(3,1))
coverageStackPlot(PTD_panel
             , alterationType=c("mutations", "copynumber")
             , var="drug"
             , grouping="gene_symbol"
             , tumor_type="brca"
             )

coverageStackPlot(PTD_panel
             , alterationType=c("mutations" , "copynumber" , "fusions")
             , var="drug"
             , grouping="gene_symbol"
             , tumor_type=c("luad", "lusc")
             )

coverageStackPlot(PTD_panel
             , alterationType=c("mutations" , "copynumber")
             , var="drug"
             , grouping="gene_symbol"
             , tumor_type=c("stad", "coadread")
             )

# reset plot settings
par(mfrow=c(1,1))

We can also visualize co-occurrence and mutual exclusivity of the drugs in different tumor types.

Overall Co-occurance

coocMutexPlot(PTD_panel
             , alterationType=c("mutations" , "copynumber")
             , var="drug"
             # , grouping="group"
             , tumor.freqs = c(brca=brca_freq
                               , luad=luad_freq
                               , lusc=lusc_freq
                               , stad = stad_freq
                               , coadread = coadread_freq
                               )
             )

Breast cancer co-occurence plot

coocMutexPlot(PTD_panel
             , alterationType=c("mutations" , "copynumber")
             #, grouping="alteration_id"
             , var="drug"
             , tumor_types= "brca"
             )

Lung cancer co-occurence plot

coocMutexPlot(PTD_panel   
             , alterationType=c("mutations" , "copynumber", "fusions" )
             #, grouping="alteration_id"
             , var="drug"
             , tumor_types= c("luad", "lusc")
             )

## Warning in dataExtractor(object = object, alterationType =
## alterationType, : The following tumor types have no alteration to display:
## stad, coadread

Gastric cancer co-occurence plot

coocMutexPlot(PTD_panel
             , alterationType=c("mutations" , "copynumber")
             #, grouping="alteration_id"
             , var="drug"
             , tumor_types= c("egc", "stad", "coadread")
             )

How many patients are eligible to more than one treatment?

Some patients are eligible to more than one therapy group. We will take advantage of this information in a later stage of the analysis as it will allows us to to optimize the patients allocation and improve the changes of the clinical trial of reaching statistical power.

# Prepare dataframe on which to run the analysis
# ------------------------------------------------------------------------------
# Extract data
data <- dataExtractor(PTD_panel , alterationType = c("mutations" , "copynumber" , "fusions"))

## Warning in dataExtractor(PTD_panel, alterationType = c("mutations",
## "copynumber", : The following tumor types have no alteration to display:
## stad, coadread

# Extract case_id and drug columns
dataDrug <- unique(data$data[ , c("case_id" , "drug")])
# calculate number of drug per patient, taking into account ALL samples
drugsPerPatient <- table( factor(dataDrug$case_id , levels = data$Samples$all_tumors) )

# Calculate the frequenecy of patients eligible to at least X number od drugs
alPatDrug <- sapply(1:max(drugsPerPatient), function(x){
          sum(drugsPerPatient >= x) / length(drugsPerPatient)
         }
       ) %>% 
  data.frame(drugs=.) 

# Plot
alPatDrug %>%
  ggplot(aes(x=1:length(drugs), y=drugs)) + 
  geom_bar(stat="identity") +
  labs(title="Patients eligible to AT LEAST 1, 2, 3 or 4 ... drugs") +
  ylab(label = "Eligibility frequency") +
  xlab(label = "Number of drugs the patient is eligible for") +
   theme(legend.position="none") + 
  geom_text(aes(label=stringr::str_c(round(drugs, digits=4)*100, " %"), y=drugs+0.05),  position = position_dodge(width=0.9))

Power analysis

We can now simulate a clinical trial based on the drug-gene targeting design, and estimate the number of patients that needs to be recruited.

• We will look at 4 power levels: 0.6, 0.7, 0.8, 0.9
• and we will consider different Hazard ratios

Overall power analysis for the full study

# Generate plot
survPowerSampleSize(PTD_panel
  , alterationType = c("mutations" , "copynumber")
  , HR = c(0.5, 0.6, 0.7, 0.8)
  , power = seq(0.6 , 0.9 , 0.1)
  , alpha = 0.05
  , case.fraction = 0.5
  , tumor.weights = c(brca=brca_tot
                    , luad = luad_tot
                    , lusc = lusc_tot
                    , stad = stad_tot
                    , coadread = coadread_tot
                    )
  )

Preview the table with the data from the previous plot.

# Generate Table
survPowerSampleSize(PTD_panel
  , HR = c(0.5, 0.6, 0.7, 0.8)
  , power = seq(0.6 , 0.9 , 0.1)
  , alpha = 0.05
  , case.fraction = 0.5
  , tumor.weights = c(brca=brca_tot
                    , luad = luad_tot
                    , lusc = lusc_tot
                    , stad = stad_tot
                    , coadread = coadread_tot
                    )
  , noPlot =TRUE
  ) %>% kable()

Var	ScreeningSampleSize	EligibleSampleSize	Beta	Power	HazardRatio
Full Design	93	59	0.4	0.6	0.5
Full Design	116	74	0.3	0.7	0.5
Full Design	149	95	0.2	0.8	0.5
Full Design	197	126	0.1	0.9	0.5
Full Design	163	104	0.4	0.6	0.6
Full Design	203	130	0.3	0.7	0.6
Full Design	259	166	0.2	0.8	0.6
Full Design	347	222	0.1	0.9	0.6
Full Design	319	204	0.4	0.6	0.7
Full Design	401	257	0.3	0.7	0.7
Full Design	511	327	0.2	0.8	0.7
Full Design	682	437	0.1	0.9	0.7
Full Design	788	505	0.4	0.6	0.8
Full Design	992	636	0.3	0.7	0.8
Full Design	1262	809	0.2	0.8	0.8
Full Design	1688	1082	0.1	0.9	0.8

Power analysis per group treatment

The power analysis can also be performed grouping by tumor types.

survPowerSampleSize(PTD_panel
  ,var = "tumor_type"
  , HR = c(0.5, 0.6, 0.7)
  , power = seq(0.6 , 0.9 , 0.1)
  , alpha = 0.05
  , case.fraction = 0.5
  , tumor.weights = c(brca=brca_tot
                    , luad = luad_tot
                    , lusc = lusc_tot
                    , stad = stad_tot
                    , coadread = coadread_tot
                    )
  )

PATIENT ALLOCATION OPTIMIZATION

In this section we will investigate how the patients could be allocated to the treatment groups. It is particularly interesting to notice that some patients are eligible to more than one therapy group, allowing to run optimization algorithms in order to better balance the patients allocation distribution in the different treatment arms reducing the number samples to screen.

Number of patients required for the overall study to reach the required power (using PFS)

In this simulation we are using Progression Free Survival as the primary endpoint.

drugs <- c("HER2 inhibitor"
  , "PARP inhibitor"
  , "AKT inhibitor"
  , "EGFR inhibitor"
  , "FGFR inhibitor"
  , "ALK inhibitor"
  , "BRAF inhibitor"
  , "MET inhibitor"
  , "NOTCH inhibitor"
  , "Immune Checkpoint Inhibitor")

# Show plot
survPowerSampleSize(PTD_panel
  , HR = c(0.5, 0.6, 0.7, 0.8)
  , power = seq(0.6 , 0.9 , 0.1)
  , alpha = 0.05
  , case.fraction = 0.5
  , priority.trial = drugs 
  , priority.trial.order="optimal"
  , var = "drug"
  , tumor.weights = c(brca=brca_tot
                    , luad = luad_tot
                    , lusc = lusc_tot
                    , stad = stad_tot
                    , coadread = coadread_tot
                    )
  )

Show in a table the Power 0.8 and HR = 0.6

# Show it in a table
survPowerSampleSize(PTD_panel
  , HR = 0.6
  , power = 0.8
  , alpha = 0.05
  , case.fraction = 0.5
  , priority.trial = drugs 
  , priority.trial.order = "optimal"
  , var = "drug"
  , noPlot = TRUE
  , tumor.weights = c(brca = brca_tot
                    , luad = luad_tot
                    , lusc = lusc_tot
                    , stad = stad_tot
                    , coadread = coadread_tot
                    )
  )$`HR:0.6 | HR0:1 | Power:0.8`$Summary %>%
  data.frame() %>%
  tibble::rownames_to_column("criteria") %>%
  tidyr::gather(drug, patients, 2:ncol(.)) %>%
  tidyr::spread(criteria, patients) %>%
  dplyr::select(drug, Screened, Eligible, Not.Eligible) %>%
  knitr::kable()

## Warning in dataExtractor(object = object, alterationType =
## alterationType, : No data for expression. They will be removed from
## alterationType

## Warning in dataExtractor(object = object, alterationType =
## alterationType, : The following tumor types have no alteration to display:
## stad, coadread

## 
## Minimum Eligible Sample Size to reach 80% of power with an HR equal to 0.6 and a HR0 equal to 1 for each drug is equal to: 166

drug	Screened	Eligible	Not.Eligible
AKT.inhibitor	0	3934	0
ALK.inhibitor	0	491	0
BRAF.inhibitor	0	397	0
EGFR.inhibitor	0	336	0
FGFR.inhibitor	0	1538	0
HER2.inhibitor	0	192	0
Immune.Checkpoint.Inhibitor	17297	166	7525
MET.inhibitor	0	286	0
NOTCH.inhibitor	0	932	0
PARP.inhibitor	0	1505	0
Total	17297	9773	7525

extra

Number of patients to be recruited to reach power

propPowerSampleSize(PTD_panel
  , noPlot=FALSE
  , alterationType = c("mutations", "copynumber")
  , pCase = c(0.2, 0.3 , 0.4)
  , pControl = rep(0.1, 3)
  , power = c(0.7, 0.8, 0.9)
  , alpha= 0.1
  , side = 1
  ) 

propPowerSampleSize(PTD_panel
  , var = "drug"
  , alterationType = c("mutations", "copynumber")
  , priority.trial = drugs 
  , priority.trial.order="optimal"
  , noPlot=TRUE
  , pCase = 0.3
  , pControl = 0.1
  , side =1
  , alpha= 0.1
  , power =0.8
  , case.fraction = 0.6
  )

## $`pCase:0.3 | pControl:0.1 | Power:0.8`
## $`pCase:0.3 | pControl:0.1 | Power:0.8`$Summary
##              Immune Checkpoint Inhibitor MET inhibitor HER2 inhibitor
## Screened                            4636             0              0
## Eligible                              65            70             84
## Not.Eligible                        2097             0              0
##              EGFR inhibitor BRAF inhibitor ALK inhibitor NOTCH inhibitor
## Screened                  0              0             0               0
## Eligible                 86            119           105             254
## Not.Eligible              0              0             0               0
##              PARP inhibitor FGFR inhibitor AKT inhibitor Total
## Screened                  0              0             0  4636
## Eligible                409            351          1001  2540
## Not.Eligible              0              0             0  2097
## 
## $`pCase:0.3 | pControl:0.1 | Power:0.8`$Screening.scheme
##                             Immune Checkpoint Inhibitor MET inhibitor
## Immune Checkpoint Inhibitor                        4636          4571
## MET inhibitor                                         0             0
## HER2 inhibitor                                        0             0
## EGFR inhibitor                                        0             0
## BRAF inhibitor                                        0             0
## ALK inhibitor                                         0             0
## NOTCH inhibitor                                       0             0
## PARP inhibitor                                        0             0
## FGFR inhibitor                                        0             0
## AKT inhibitor                                         0             0
##                             HER2 inhibitor EGFR inhibitor BRAF inhibitor
## Immune Checkpoint Inhibitor           4502           4418           4332
## MET inhibitor                            0              0              0
## HER2 inhibitor                           0              0              0
## EGFR inhibitor                           0              0              0
## BRAF inhibitor                           0              0              0
## ALK inhibitor                            0              0              0
## NOTCH inhibitor                          0              0              0
## PARP inhibitor                           0              0              0
## FGFR inhibitor                           0              0              0
## AKT inhibitor                            0              0              0
##                             ALK inhibitor NOTCH inhibitor PARP inhibitor
## Immune Checkpoint Inhibitor          4214            4110           3856
## MET inhibitor                           0               0              0
## HER2 inhibitor                          0               0              0
## EGFR inhibitor                          0               0              0
## BRAF inhibitor                          0               0              0
## ALK inhibitor                           0               0              0
## NOTCH inhibitor                         0               0              0
## PARP inhibitor                          0               0              0
## FGFR inhibitor                          0               0              0
## AKT inhibitor                           0               0              0
##                             FGFR inhibitor AKT inhibitor
## Immune Checkpoint Inhibitor           3448          3097
## MET inhibitor                            0             0
## HER2 inhibitor                           0             0
## EGFR inhibitor                           0             0
## BRAF inhibitor                           0             0
## ALK inhibitor                            0             0
## NOTCH inhibitor                          0             0
## PARP inhibitor                           0             0
## FGFR inhibitor                           0             0
## AKT inhibitor                            0             0
## 
## $`pCase:0.3 | pControl:0.1 | Power:0.8`$Allocation.scheme
##                             Immune Checkpoint Inhibitor MET inhibitor
## Immune Checkpoint Inhibitor                          65            70
## MET inhibitor                                         0             0
## HER2 inhibitor                                        0             0
## EGFR inhibitor                                        0             0
## BRAF inhibitor                                        0             0
## ALK inhibitor                                         0             0
## NOTCH inhibitor                                       0             0
## PARP inhibitor                                        0             0
## FGFR inhibitor                                        0             0
## AKT inhibitor                                         0             0
## Total                                                65            70
##                             HER2 inhibitor EGFR inhibitor BRAF inhibitor
## Immune Checkpoint Inhibitor             84             86            119
## MET inhibitor                            0              0              0
## HER2 inhibitor                           0              0              0
## EGFR inhibitor                           0              0              0
## BRAF inhibitor                           0              0              0
## ALK inhibitor                            0              0              0
## NOTCH inhibitor                          0              0              0
## PARP inhibitor                           0              0              0
## FGFR inhibitor                           0              0              0
## AKT inhibitor                            0              0              0
## Total                                   84             86            119
##                             ALK inhibitor NOTCH inhibitor PARP inhibitor
## Immune Checkpoint Inhibitor           105             254            409
## MET inhibitor                           0               0              0
## HER2 inhibitor                          0               0              0
## EGFR inhibitor                          0               0              0
## BRAF inhibitor                          0               0              0
## ALK inhibitor                           0               0              0
## NOTCH inhibitor                         0               0              0
## PARP inhibitor                          0               0              0
## FGFR inhibitor                          0               0              0
## AKT inhibitor                           0               0              0
## Total                                 105             254            409
##                             FGFR inhibitor AKT inhibitor
## Immune Checkpoint Inhibitor            351          1001
## MET inhibitor                            0             0
## HER2 inhibitor                           0             0
## EGFR inhibitor                           0             0
## BRAF inhibitor                           0             0
## ALK inhibitor                            0             0
## NOTCH inhibitor                          0             0
## PARP inhibitor                           0             0
## FGFR inhibitor                           0             0
## AKT inhibitor                            0             0
## Total                                  351          1001
## 
## $`pCase:0.3 | pControl:0.1 | Power:0.8`$Probability.scheme
##                             Immune Checkpoint Inhibitor MET inhibitor
## Immune Checkpoint Inhibitor                  0.01402103    0.01523616
## MET inhibitor                                0.00000000    0.00000000
## HER2 inhibitor                               0.00000000    0.00000000
## EGFR inhibitor                               0.00000000    0.00000000
## BRAF inhibitor                               0.00000000    0.00000000
## ALK inhibitor                                0.00000000    0.00000000
## NOTCH inhibitor                              0.00000000    0.00000000
## PARP inhibitor                               0.00000000    0.00000000
## FGFR inhibitor                               0.00000000    0.00000000
## AKT inhibitor                                0.00000000    0.00000000
##                             HER2 inhibitor EGFR inhibitor BRAF inhibitor
## Immune Checkpoint Inhibitor     0.01856627     0.01944298     0.02733119
## MET inhibitor                   0.00000000     0.00000000     0.00000000
## HER2 inhibitor                  0.00000000     0.00000000     0.00000000
## EGFR inhibitor                  0.00000000     0.00000000     0.00000000
## BRAF inhibitor                  0.00000000     0.00000000     0.00000000
## ALK inhibitor                   0.00000000     0.00000000     0.00000000
## NOTCH inhibitor                 0.00000000     0.00000000     0.00000000
## PARP inhibitor                  0.00000000     0.00000000     0.00000000
## FGFR inhibitor                  0.00000000     0.00000000     0.00000000
## AKT inhibitor                   0.00000000     0.00000000     0.00000000
##                             ALK inhibitor NOTCH inhibitor PARP inhibitor
## Immune Checkpoint Inhibitor    0.02479339      0.06158192      0.1059603
## MET inhibitor                  0.00000000      0.00000000      0.0000000
## HER2 inhibitor                 0.00000000      0.00000000      0.0000000
## EGFR inhibitor                 0.00000000      0.00000000      0.0000000
## BRAF inhibitor                 0.00000000      0.00000000      0.0000000
## ALK inhibitor                  0.00000000      0.00000000      0.0000000
## NOTCH inhibitor                0.00000000      0.00000000      0.0000000
## PARP inhibitor                 0.00000000      0.00000000      0.0000000
## FGFR inhibitor                 0.00000000      0.00000000      0.0000000
## AKT inhibitor                  0.00000000      0.00000000      0.0000000
##                             FGFR inhibitor AKT inhibitor
## Immune Checkpoint Inhibitor      0.1016835     0.3230885
## MET inhibitor                    0.0000000     0.0000000
## HER2 inhibitor                   0.0000000     0.0000000
## EGFR inhibitor                   0.0000000     0.0000000
## BRAF inhibitor                   0.0000000     0.0000000
## ALK inhibitor                    0.0000000     0.0000000
## NOTCH inhibitor                  0.0000000     0.0000000
## PARP inhibitor                   0.0000000     0.0000000
## FGFR inhibitor                   0.0000000     0.0000000
## AKT inhibitor                    0.0000000     0.0000000
## 
## $`pCase:0.3 | pControl:0.1 | Power:0.8`$Base.probabilities
## Immune Checkpoint Inhibitor               MET inhibitor 
##                  0.01402103                  0.01602404 
##              HER2 inhibitor              EGFR inhibitor 
##                  0.02003005                  0.02153230 
##              BRAF inhibitor               ALK inhibitor 
##                  0.03054582                  0.03254882 
##             NOTCH inhibitor              PARP inhibitor 
##                  0.06659990                  0.10916375 
##              FGFR inhibitor               AKT inhibitor 
##                  0.11417126                  0.31347021

surv2 <- survPowerSampleSize(PTD_panel
   , alterationType = c("mutations", "copynumber")
   , HR = 0.66
   , power = 0.8
   , alpha = 0.05
   , case.fraction = 0.5
   , ber=0.85
   , fu=1.883
   #, priority.trial.order="as.is"
   , priority.trial.order="optimal"
   , var = "drug" 
   , noPlot=TRUE)
kable(surv2 , row.names=FALSE)

Var	ScreeningSampleSize	EligibleSampleSize	Beta	Power	HazardRatio
AKT inhibitor	801	251	0.2	0.8	0.66
ALK inhibitor	7712	251	0.2	0.8	0.66
BRAF inhibitor	8218	251	0.2	0.8	0.66
EGFR inhibitor	11657	251	0.2	0.8	0.66
FGFR inhibitor	2199	251	0.2	0.8	0.66
HER2 inhibitor	12532	251	0.2	0.8	0.66
Immune Checkpoint Inhibitor	17902	251	0.2	0.8	0.66
MET inhibitor	15664	251	0.2	0.8	0.66
no_drug	1700	251	0.2	0.8	0.66
NOTCH inhibitor	3769	251	0.2	0.8	0.66
PARP inhibitor	2300	251	0.2	0.8	0.66
Full Design	398	251	0.2	0.8	0.66

For the Abstract

PFS 2 samples calculations

Perform a 1-side 2-sample sample size calculation for PFS study using:

• fu 1.883
• case.fraction = 0.5
• ber = 0.85
• HR = 0.625
• power = 0.8

NOT OPTIMAL

library(parallel)
set.seed(111)
# NOT optimal
surv2 <- mclapply(1:100 , function(x){
  out <- survPowerSampleSize(PTD_panel
   , alterationType = c("mutations", "copynumber")
   , HR = 0.625
   , power = 0.8
   , alpha = 0.05
   , case.fraction = 0.5
   , ber=0.85
   , fu=1.883
   , var = "drug"
   , side = 1
   , tumor.weights = c(brca=brca_tot
                               , luad=luad_tot
                               , lusc=lusc_tot
                               , stad = stad_tot
                               , coadread = coadread_tot
                               )
   , noPlot=TRUE)
  out <- out[ order(out$Var) , ]
  return(out)
} , mc.cores=detectCores())

Screenings <- sapply(surv2 , function(x) x$ScreeningSampleSize)
medianScreening <- apply(Screenings , 1 , function(x) {
  paste(median(x) 
        , paste0( "(" , quantile(x , 0.025) 
                        , "-" 
                        , quantile(x , 0.975) , ")"))
})

surv2DF <- data.frame(Var = surv2[[1]]$Var
                  , ScreeningSampleSize = medianScreening
                  , EligibleSampleSize = surv2[[1]]$EligibleSampleSize
                  , Beta = 0.2
                  , Power = 0.8
                  , HR = 0.625)
#export
#write.table(surv2, file = "~/Desktop/notoptimal_ptd-showcase.txt", quote=FALSE, sep=",")

#preview
kable(surv2DF[ match(c("Full Design" , drugs , "no_drug") , surv2DF$Var) , ] , row.names=FALSE)

Var	ScreeningSampleSize	EligibleSampleSize	Beta	Power	HR
Full Design	241 (237-245.525)	157	0.2	0.8	0.625
HER2 inhibitor	9073.5 (7667.25-11185.125)	157	0.2	0.8	0.625
PARP inhibitor	1748 (1659.425-1856.25)	157	0.2	0.8	0.625
AKT inhibitor	529 (507.475-546.575)	157	0.2	0.8	0.625
EGFR inhibitor	6327 (5495-7455)	157	0.2	0.8	0.625
FGFR inhibitor	1595 (1492.125-1725.925)	157	0.2	0.8	0.625
ALK inhibitor	4615 (4076.6-5101)	157	0.2	0.8	0.625
BRAF inhibitor	3767 (3383.425-4118.075)	157	0.2	0.8	0.625
MET inhibitor	9434 (8135.375-11654.675)	157	0.2	0.8	0.625
NOTCH inhibitor	2798 (2584.025-3046.925)	157	0.2	0.8	0.625
Immune Checkpoint Inhibitor	11845 (9963-14314.8)	157	0.2	0.8	0.625
no_drug	731 (707.475-759.575)	157	0.2	0.8	0.625

• Median Number of patients necessary with CI: 52290 (48562.95-55701.02)

PRIORITY OPTIMAL

set.seed(112)
simulSurv2 <- mclapply(1:100 , function(x) {
    survPowerSampleSize(PTD_panel
        , alterationType = c("mutations", "copynumber")
        , HR = 0.625
        , power = 0.8
        , alpha = 0.05
        , ber=0.85
        , fu=1.883
        , case.fraction = 0.5
        , collapseByGene=TRUE
        , tumor.weights = c(brca=brca_tot
                               , luad=luad_tot
                               , lusc=lusc_tot
                               , stad = stad_tot
                               , coadread = coadread_tot
                               )
        , var = "drug" 
        , noPlot=TRUE
        , side = 1
        , priority.trial=drugs
        , priority.trial.order="optimal")$`HR:0.625 | HR0:1 | Power:0.8`$Summary
} , mc.cores=detectCores())
summaryMat2 <- simulSurv2[[1]]
for(i in rownames(summaryMat2)){
    for(j in colnames(summaryMat2)){
        vec <- sapply(simulSurv2 , function(x) x[i , j])
        summaryMat2[i , j] <- paste(median(vec) , 
                paste0( "(" , quantile(vec , 0.025) 
                        , "-" 
                        , quantile(vec , 0.975) , ")"))
    }
}

#save
#write.table(summaryMat2, file = "~/Desktop/optimal_ptd-showcase.txt", quote=FALSE, sep=",")
#preview
kable(summaryMat2)

	Immune Checkpoint Inhibitor	MET inhibitor	HER2 inhibitor	EGFR inhibitor	ALK inhibitor	BRAF inhibitor	NOTCH inhibitor	PARP inhibitor	FGFR inhibitor	AKT inhibitor	Total
Screened	11845 (639.4-14125.175)	0 (0-11121.075)	0 (0-716.449999999997)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	12114 (10382.725-14125.175)
Eligible	157 (157-171.05)	178 (157-222.05)	192 (157-224.05)	276.5 (227.9-346)	316 (252.95-370.05)	414 (347.9-499.1)	559 (456.55-656)	856.5 (722.75-1007.25)	807 (680.55-978.825)	2514.5 (2142.275-2979.075)	6272.5 (5413.775-7355.125)
Not.Eligible	5707 (308.65-6866.9)	0 (0-5323.8)	0 (0-345.174999999998)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	5762 (4969.475-6866.9)

Number of patients necessary under optimal design (median value with CI): 12114 (10382.725-14125.175)

Basically the number of samples screened for Immune check inhibitor, that is targetable only in a very small number of patients, is sufficient to cover all the other drugs in most of the cases.

The percentage of samples saved from passing from a non optimal (all drugs separated) to an optimal design (all drugs evaluated together) is 76.83%

PRIORITY NON OPTIMAL

Extract drug frequencies

drugFreqs <- coverageStackPlot(PTD_panel
   , alterationType = c("mutations", "copynumber")
   , var = "drug"
   , grouping = NA
   , tumor.freqs = c(brca=brca_freq
                    , luad=luad_freq
                    , lusc=lusc_freq
                    , stad = stad_freq
                    , coadread = coadread_freq) 
        , noPlot=TRUE)$plottedTable[1, ]
# Sort from highest to lowest
drugFreqs <- sort(drugFreqs , decreasing = TRUE)
# Remove no_drug
drugFreqs <- drugFreqs[ names(drugFreqs)!="no_drug"]

Perform a priority trial starting from the drug with the highest frequency (NON OPTIMAL)

set.seed(112)
simulSurv2_NO <- mclapply(1:100 , function(x) {
    survPowerSampleSize(PTD_panel
        , alterationType = c("mutations", "copynumber")
        , HR = 0.625
        , power = 0.8
        , alpha = 0.05
        , ber=0.85
        , fu=1.883
        , case.fraction = 0.5
        , collapseByGene=TRUE
        , tumor.weights = c(brca=brca_tot
                               , luad=luad_tot
                               , lusc=lusc_tot
                               , stad = stad_tot
                               , coadread = coadread_tot
                               )
        , var = "drug" 
        , noPlot=TRUE
        , side = 1
        , priority.trial=names(drugFreqs)
        , priority.trial.order="as.is")$`HR:0.625 | HR0:1 | Power:0.8`$Summary
} , mc.cores=detectCores())
summaryMat2_NO <- simulSurv2_NO[[1]]
for(i in rownames(summaryMat2_NO)){
    for(j in colnames(summaryMat2_NO)){
        vec <- sapply(simulSurv2_NO , function(x) x[i , j])
        summaryMat2_NO[i , j] <- paste(median(vec) , 
                paste0( "(" , quantile(vec , 0.025) 
                        , "-" 
                        , quantile(vec , 0.975) , ")"))
    }
}

#save
#write.table(summaryMat2_NO, file = "~/Desktop/optimal_ptd-showcase.txt", quote=FALSE, sep=",")
#preview
kable(summaryMat2_NO)

	AKT inhibitor	FGFR inhibitor	PARP inhibitor	NOTCH inhibitor	BRAF inhibitor	ALK inhibitor	EGFR inhibitor	HER2 inhibitor	MET inhibitor	Immune Checkpoint Inhibitor	Total
Screened	525 (509.475-544.05)	1220 (1123.275-1330.325)	654.5 (480.75-821.625)	1737 (1441.475-2043.825)	1064.5 (593.125-1617.05)	2424.5 (1956.225-3230)	2343 (1348.875-3440.2)	5956.5 (4384.45-7857.55)	2891.5 (573.7-5285.8)	5939 (3597.325-8964.25)	24812 (22029.175-28888.525)
Eligible	7121 (6338.1-8300.45)	1717.5 (1508.5-2041.425)	1380.5 (1187.275-1664)	674 (567.225-822.35)	677 (561.85-794.675)	398 (323.85-489.05)	346.5 (300.475-423.525)	205 (184.475-236.05)	201 (180.95-236)	157 (157-157)	12861 (11491.975-14980.225)
Not.Eligible	252 (239.475-266)	585.5 (533.425-635.05)	312.5 (232.475-400.575)	819.5 (691.325-979.775)	505 (279.75-766.175)	1153.5 (919.025-1545.45)	1123.5 (640.05-1653)	2843 (2073.55-3813.875)	1391 (270.875-2526.6)	5859 (3543.225-8853.775)	14991.5 (12854.15-18217.875)

Number of patients necessary under non optimal priority design (median value with CI): 24812 (22029.175-28888.525)

FULL SUMMARY OF THE THREE METHODS

Create a dataframe with median and CI for the three methods

pfs2 <- data.frame( Method = c("Indipendent" , "Priority Non Optimal" , "Priority Optimal")
                    , Median = c(median(n_patients_notoptimal) 
                                 , median(n_patients_optimal_NO) 
                                 , median(n_patients_optimal))
                    , CI_low = c(quantile(n_patients_notoptimal , 0.025) 
                                 , quantile(n_patients_optimal_NO , 0.025) 
                                 , quantile(n_patients_optimal , 0.025))
                    , CI_high = c(quantile(n_patients_notoptimal , 0.975) 
                                 , quantile(n_patients_optimal_NO , 0.975) 
                                 , quantile(n_patients_optimal , 0.975)))

Plot as barplot

gpfs2 <- ggplot(pfs2, aes(x=Method, y=Median)) + 
    geom_bar(position=position_dodge(), stat="identity" 
             , fill = c("firebrick3" , "deepskyblue4" , "seagreen4") 
            , col=rep("black",3)) +
    geom_errorbar(aes(ymin=CI_low, ymax=CI_high)
                  , width=.2
                  , position=position_dodge(.9)
                  , size=0.2
                  )
gpfs2

ggsave(filename="../Figures/fig5A.svg", plot=gpfs2, device = "svg")

## Saving 10 x 8 in image

PFS 1 sample calculations

Perform a 1-side 1-sample sample size calculation for PFS study using:

• MED0 = 0.8154673 (corresponding to ber = 0.85)
• MED1 = 1.304748 (corresponding to MED0/HR)
• power = 0.8
• fu = 1.883

NOT OPTIMAL

set.seed(111)
# NOT optimal
surv3 <- mclapply(1:100 , function(x){
  out <- survPowerSampleSize1Arm(PTD_panel
   , alterationType = c("mutations", "copynumber")
   , MED0 = 0.8154673
   , MED1 = 1.304748
   , power = 0.8
   , alpha = 0.05
   , fu = 1.883
   , var = "drug"
   , side = 1
   , tumor.weights = c(brca=brca_tot
                               , luad=luad_tot
                               , lusc=lusc_tot
                               , stad = stad_tot
                               , coadread = coadread_tot
                               )
   , noPlot=TRUE)
  out <- out[ order(out$Var) , ]
  return(out)
} , mc.cores=detectCores())

Screenings <- sapply(surv3 , function(x) x$ScreeningSampleSize)
medianScreening <- apply(Screenings , 1 , function(x) {
  paste(median(x) 
        , paste0( "(" , quantile(x , 0.025) 
                        , "-" 
                        , quantile(x , 0.975) , ")"))
})

surv3DF <- data.frame(Var = surv3[[1]]$Var
                  , ScreeningSampleSize = medianScreening
                  , EligibleSampleSize = surv3[[1]]$EligibleSampleSize
                  , Beta = 0.2
                  , Power = 0.8
                  , MED1 = 1.304748)
#export
#write.table(surv3, file = "~/Desktop/notoptimal_ptd-showcase.txt", quote=FALSE, sep=",")

#preview
kable(surv3DF[ match(c("Full Design" , drugs , "no_drug") , surv3DF$Var) , ] , row.names=FALSE)

Var	ScreeningSampleSize	EligibleSampleSize	Beta	Power	MED1
Full Design	62 (61-63)	40	0.2	0.8	1.304748
HER2 inhibitor	2312 (1983-2718.35)	40	0.2	0.8	1.304748
PARP inhibitor	454.5 (413.95-484.2)	40	0.2	0.8	1.304748
AKT inhibitor	135 (130-139)	40	0.2	0.8	1.304748
EGFR inhibitor	1617 (1425.8-1887)	40	0.2	0.8	1.304748
FGFR inhibitor	407 (379.475-438.05)	40	0.2	0.8	1.304748
ALK inhibitor	1184 (1058.9-1332)	40	0.2	0.8	1.304748
BRAF inhibitor	953 (849.375-1055.3)	40	0.2	0.8	1.304748
MET inhibitor	2492 (2097.6-2937.8)	40	0.2	0.8	1.304748
NOTCH inhibitor	718 (673.475-778.05)	40	0.2	0.8	1.304748
Immune Checkpoint Inhibitor	3018 (2539-3647.725)	40	0.2	0.8	1.304748
no_drug	187 (179-192)	40	0.2	0.8	1.304748

• Median Number of patients necessary with CI: 13307 (12390.2-14261.5)

OPTIMAL

set.seed(112)
simulSurv3 <- mclapply(1:100 , function(x) {
    survPowerSampleSize1Arm(PTD_panel
    , alterationType = c("mutations", "copynumber")
    , MED0 = 0.8154673
    , MED1 = 1.304748
    , power = 0.8
    , alpha = 0.05
    , fu = 1.883
    , var = "drug"
    , side = 1
    , tumor.weights = c(brca=brca_tot
                               , luad=luad_tot
                               , lusc=lusc_tot
                               , stad = stad_tot
                               , coadread = coadread_tot
                               )
    , noPlot=TRUE
    , priority.trial=drugs
    , priority.trial.order="optimal")$`MED1:1.304748 | MED0:0.8154673 | Power:0.8`$Summary
} , mc.cores=detectCores())
summaryMat3 <- simulSurv3[[1]]
for(i in rownames(summaryMat3)){
    for(j in colnames(summaryMat3)){
        vec <- sapply(simulSurv3 , function(x) x[i , j])
        summaryMat3[i , j] <- paste(median(vec) , 
                paste0( "(" , quantile(vec , 0.025) 
                        , "-" 
                        , quantile(vec , 0.975) , ")"))
    }
}

#save
#write.table(summaryMat3, file = "~/Desktop/optimal_ptd-showcase.txt", quote=FALSE, sep=",")
#preview
kable(summaryMat3)

	Immune Checkpoint Inhibitor	HER2 inhibitor	MET inhibitor	EGFR inhibitor	ALK inhibitor	BRAF inhibitor	NOTCH inhibitor	PARP inhibitor	FGFR inhibitor	AKT inhibitor	Total
Screened	3018 (274-3574)	0 (0-200.825)	0 (0-1497.34999999998)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	3018 (2649.35-3577.15)
Eligible	40 (40-43.525)	47.5 (40-63.05)	46 (40-58.525)	69 (56.95-84)	79 (64-97.05)	104.5 (83-123.525)	139.5 (119.95-174.575)	214 (181.9-266)	204.5 (176.95-246.525)	637 (546.375-772.15)	1580.5 (1374.65-1879.775)
Not.Eligible	1437.5 (130.475-1715.15)	0 (0-96.1249999999998)	0 (0-718.024999999991)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	1447 (1266.325-1715.15)

Number of patients necessary under optimal design (median value with CI): 3018 (2649.35-3577.15)

Basically the number of samples screened for Immune check inhibitor, that is targetable only in a very small number of patients, is sufficient to cover all the other drugs in most of the cases.

The percentage of samples saved from passing from a non optimal (all drugs separated) to an optimal design (all drugs evaluated together) is 77.32%

PRIORITY NON OPTIMAL

Perform a priority trial starting from the drug with the highest frequency (NON OPTIMAL)

set.seed(112)
simulSurv3_NO <- mclapply(1:100 , function(x) {
    survPowerSampleSize1Arm(PTD_panel
    , alterationType = c("mutations", "copynumber")
    , MED0 = 0.8154673
    , MED1 = 1.304748
    , power = 0.8
    , alpha = 0.05
    , fu = 1.883
    , var = "drug"
    , side = 1
    , tumor.weights = c(brca=brca_tot
                        , luad=luad_tot
                        , lusc=lusc_tot
                        , stad = stad_tot
                        , coadread = coadread_tot)
    , noPlot=TRUE
    , priority.trial=names(drugFreqs)
    , priority.trial.order="as.is")$`MED1:1.304748 | MED0:0.8154673 | Power:0.8`$Summary
} , mc.cores=detectCores())
summaryMat3_NO <- simulSurv3_NO[[1]]
for(i in rownames(summaryMat3_NO)){
    for(j in colnames(summaryMat3_NO)){
        vec <- sapply(simulSurv3_NO , function(x) x[i , j])
        summaryMat3_NO[i , j] <- paste(median(vec) , 
                paste0( "(" , quantile(vec , 0.025) 
                        , "-" 
                        , quantile(vec , 0.975) , ")"))
    }
}

#save
#write.table(summaryMat3, file = "~/Desktop/optimal_ptd-showcase.txt", quote=FALSE, sep=",")
#preview
kable(summaryMat3_NO)

	AKT inhibitor	FGFR inhibitor	PARP inhibitor	NOTCH inhibitor	BRAF inhibitor	ALK inhibitor	EGFR inhibitor	HER2 inhibitor	MET inhibitor	Immune Checkpoint Inhibitor	Total
Screened	135 (131-139)	307.5 (282.425-343.15)	172 (126.95-220)	448.5 (381.475-523.2)	259.5 (132.9-398.2)	646.5 (491.225-840.575)	560.5 (269.9-819.8)	1457.5 (1098.1-1954.2)	674 (83.275-1328.275)	1670 (1039.35-2362.975)	6339.5 (5668.8-7159.95)
Eligible	1822 (1603.525-2098.15)	444 (388.375-530.6)	349.5 (299.975-404.575)	171 (139.475-201.1)	175 (147.475-210.15)	99 (83-121)	92 (77.475-111.525)	54 (47.475-60.525)	53 (47-60.575)	40 (40-40)	3290 (2948.2-3750.575)
Not.Eligible	65.5 (63-69)	147.5 (134.95-166.1)	82 (61.475-105.575)	215.5 (182.425-252.05)	125 (63.95-189.625)	309 (230.9-406)	264.5 (129.9-393.775)	695 (525.225-944.675)	328 (39.975-645.65)	1648.5 (1023.4-2336.5)	3859 (3221.55-4638.55)

Number of patients necessary under non optimal priority design (median value with CI): 6339.5 (5668.8-7159.95)

FULL SUMMARY OF THE THREE METHODS

Create a dataframe with median and CI for the three methods

pfs1 <- data.frame( Method = c("Indipendent" , "Priority Non Optimal" , "Priority Optimal")
                    , Median = c(median(n_patients_notoptimal) 
                                 , median(n_patients_optimal_NO) 
                                 , median(n_patients_optimal))
                    , CI_low = c(quantile(n_patients_notoptimal , 0.025) 
                                 , quantile(n_patients_optimal_NO , 0.025) 
                                 , quantile(n_patients_optimal , 0.025))
                    , CI_high = c(quantile(n_patients_notoptimal , 0.975) 
                                 , quantile(n_patients_optimal_NO , 0.975) 
                                 , quantile(n_patients_optimal , 0.975)))
kable(pfs1)

Method	Median	CI_low	CI_high
Indipendent	13307.0	12390.20	14261.50
Priority Non Optimal	6339.5	5668.80	7159.95
Priority Optimal	3018.0	2649.35	3577.15

Plot as barplot

gpfs1 <- ggplot(pfs1, aes(x=Method, y=Median)) + 
    geom_bar(position=position_dodge(), stat="identity" 
             , fill = c("firebrick3" , "deepskyblue4" , "seagreen4") 
            , col=rep("black",3)) +
    geom_errorbar(aes(ymin=CI_low, ymax=CI_high)
                  , width=.2
                  , position=position_dodge(.9)
                  , size=0.2
                  )
gpfs1

ggsave(filename="../Figures/fig5B.svg", plot=gpfs1, device = "svg")

## Saving 10 x 8 in image

Proportion trial 2-sample

Perform a 1-side 2-sample sample size calculation for difference in proportion with:

• pCase = 0.4
• pControl = 0.1
• case.fraction = 0.5
• power = 0.8

NOT OPTIMAL

set.seed(111)
# NOT optimal
prop2 <- mclapply(1:100 , function(x){
  out <- propPowerSampleSize(PTD_panel
   , alterationType = c("mutations", "copynumber")
   , pCase = 0.4
   , pControl = 0.1
   , power = 0.8
   , case.fraction = 0.5
   , type = "chisquare"
   , var = "drug"
   , side = 1
   , tumor.weights = c(brca=brca_tot
                               , luad=luad_tot
                               , lusc=lusc_tot
                               , stad = stad_tot
                               , coadread = coadread_tot
                               )
   , noPlot=TRUE)
  out <- out[ order(out$Var) , ]
  return(out)
} , mc.cores=detectCores())

Screenings <- sapply(prop2 , function(x) x$ScreeningSampleSize)
medianScreening <- apply(Screenings , 1 , function(x) {
  paste(median(x) 
        , paste0( "(" , quantile(x , 0.025) 
                        , "-" 
                        , quantile(x , 0.975) , ")"))
})

prop2DF <- data.frame(Var = prop2[[1]]$Var
                  , ScreeningSampleSize = medianScreening
                  , EligibleSampleSize = prop2[[1]]$EligibleSampleSize
                  , Beta = 0.2
                  , Power = 0.8
                  , Proportion.In.Case.Control = "0.4 - 0.1")
#export
#write.table(surv2, file = "~/Desktop/notoptimal_ptd-showcase.txt", quote=FALSE, sep=",")

#preview
kable(prop2DF[ match(c("Full Design" , drugs , "no_drug") , prop2DF$Var) , ] , row.names=FALSE)

Var	ScreeningSampleSize	EligibleSampleSize	Beta	Power	Proportion.In.Case.Control
Full Design	71 (70-72)	46	0.2	0.8	0.4 - 0.1
HER2 inhibitor	2670 (2246.4-3140.275)	46	0.2	0.8	0.4 - 0.1
PARP inhibitor	517 (476.9-557.35)	46	0.2	0.8	0.4 - 0.1
AKT inhibitor	155 (149-159.525)	46	0.2	0.8	0.4 - 0.1
EGFR inhibitor	1887.5 (1675.275-2200)	46	0.2	0.8	0.4 - 0.1
FGFR inhibitor	467.5 (442.475-501.525)	46	0.2	0.8	0.4 - 0.1
ALK inhibitor	1355.5 (1163.075-1521.075)	46	0.2	0.8	0.4 - 0.1
BRAF inhibitor	1098 (980.675-1237.625)	46	0.2	0.8	0.4 - 0.1
MET inhibitor	2827 (2339.175-3417.8)	46	0.2	0.8	0.4 - 0.1
NOTCH inhibitor	823.5 (751.95-919.775)	46	0.2	0.8	0.4 - 0.1
Immune Checkpoint Inhibitor	3510 (2932.825-4194.4)	46	0.2	0.8	0.4 - 0.1
no_drug	214 (206.425-223)	46	0.2	0.8	0.4 - 0.1

• Median Number of patients necessary with CI: 15294 (14306.45-16671.9)

OPTIMAL

set.seed(112)
simulprop2 <- mclapply(1:100 , function(x) {
    propPowerSampleSize(PTD_panel
   , alterationType = c("mutations", "copynumber")
   , pCase = 0.4
   , pControl = 0.1
   , power = 0.8
   , case.fraction = 0.5
   , var = "drug"
   , type = "chisquare"
   , tumor.weights = c(brca=brca_tot
                               , luad=luad_tot
                               , lusc=lusc_tot
                               , stad = stad_tot
                               , coadread = coadread_tot
                               ) 
        , noPlot=TRUE
        , priority.trial=drugs
        , priority.trial.order="optimal")$`pCase:0.4 | pControl:0.1 | Power:0.8`$Summary
} , mc.cores=detectCores())
prop2summaryMat2 <- simulprop2[[1]]
for(i in rownames(prop2summaryMat2)){
    for(j in colnames(prop2summaryMat2)){
        vec <- sapply(simulprop2 , function(x) x[i , j])
        prop2summaryMat2[i , j] <- paste(median(vec) , 
                paste0( "(" , quantile(vec , 0.025) 
                        , "-" 
                        , quantile(vec , 0.975) , ")"))
    }
}

#save
#write.table(summaryMat2, file = "~/Desktop/optimal_ptd-showcase.txt", quote=FALSE, sep=",")
#preview
kable(prop2summaryMat2)

	Immune Checkpoint Inhibitor	MET inhibitor	HER2 inhibitor	EGFR inhibitor	ALK inhibitor	BRAF inhibitor	NOTCH inhibitor	PARP inhibitor	FGFR inhibitor	AKT inhibitor	Total
Screened	4376 (250.525-5322)	0 (0-3901.95)	0 (0-2372.32499999997)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	4476 (3918.525-5322)
Eligible	58 (58-63)	66 (58-85)	68 (58-86.525)	104 (88.475-128.1)	117.5 (95.95-147.05)	152.5 (130.475-193.1)	207 (172.425-256.625)	319 (277.95-366)	295.5 (251.425-353.05)	945.5 (821.375-1125.775)	2326 (2058.8-2725.25)
Not.Eligible	2100 (118.1-2598.25)	0 (0-1882.425)	0 (0-1136.07499999999)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	2148.5 (1877.275-2598.25)

Number of patients necessary under optimal design (median value): 4476 (3918.525-5322)

Basically the number of samples screened for Immune check inhibitor, that is targetable only in a very small number of patients, is sufficient to cover all the other drugs in most of the cases.

The percentage of samples saved from passing from a non optimal (all drugs separated) to an optimal design (all drugs evaluated together) is 70.73%

PRIORITY NON OPTIMAL

Perform a priority trial starting from the drug with the highest frequency (NON OPTIMAL)

set.seed(112)
simulprop2_NO <- mclapply(1:100 , function(x) {
    propPowerSampleSize(PTD_panel
   , alterationType = c("mutations", "copynumber")
   , pCase = 0.4
   , pControl = 0.1
   , power = 0.8
   , case.fraction = 0.5
   , var = "drug"
   , type = "chisquare"
   , tumor.weights = c(brca=brca_tot
                               , luad=luad_tot
                               , lusc=lusc_tot
                               , stad = stad_tot
                               , coadread = coadread_tot
                               ) 
        , noPlot=TRUE
        , priority.trial=names(drugFreqs)
        , priority.trial.order="as.is")$`pCase:0.4 | pControl:0.1 | Power:0.8`$Summary
} , mc.cores=detectCores())
prop2summaryMat2_NO <- simulprop2_NO[[1]]
for(i in rownames(prop2summaryMat2_NO)){
    for(j in colnames(prop2summaryMat2_NO)){
        vec <- sapply(simulprop2_NO , function(x) x[i , j])
        prop2summaryMat2_NO[i , j] <- paste(median(vec) , 
                paste0( "(" , quantile(vec , 0.025) 
                        , "-" 
                        , quantile(vec , 0.975) , ")"))
    }
}
#save
#write.table(summaryMat3, file = "~/Desktop/optimal_ptd-showcase.txt", quote=FALSE, sep=",")
#preview
kable(prop2summaryMat2_NO)

	AKT inhibitor	FGFR inhibitor	PARP inhibitor	NOTCH inhibitor	BRAF inhibitor	ALK inhibitor	EGFR inhibitor	HER2 inhibitor	MET inhibitor	Immune Checkpoint Inhibitor	Total
Screened	195 (188-201.525)	450 (411.475-489.1)	244 (187.85-297.575)	642 (529.8-760.4)	384.5 (230.95-579.2)	930 (711.6-1131.1)	861.5 (372.3-1416.95)	2200.5 (1741.575-2841.775)	1030 (111.125-1707.175)	2311.5 (1321.625-3622.475)	9402 (7803.425-10717.6)
Eligible	2703.5 (2242.025-3094.1)	655.5 (547.8-775.625)	520 (417.275-618.1)	249 (214-299)	253.5 (214-305)	148 (121-169.575)	132 (107.425-156)	76 (69-86.525)	74.5 (67-85)	58 (58-59)	4880.5 (4149.8-5555.9)
Not.Eligible	94 (89-98)	215.5 (195-237)	118 (89.9-144.05)	309 (252.425-362)	184 (110.475-281.525)	446 (341.375-550.3)	412.5 (174.375-682.025)	1056 (839.45-1395.6)	489.5 (54.325-832.15)	2280.5 (1300.725-3582.625)	5643 (4441.425-6888.875)

Number of patients necessary under non optimal priority design (median value with CI): 6339.5 (5668.8-7159.95)

FULL SUMMARY OF THE THREE METHODS

Create a dataframe with median and CI for the three methods

propd <- data.frame( Method = c("Indipendent" , "Priority Non Optimal" , "Priority Optimal")
                    , Median = c(median(n_patients_notoptimal) 
                                 , median(n_patients_optimal_NO) 
                                 , median(n_patients_optimal))
                    , CI_low = c(quantile(n_patients_notoptimal , 0.025) 
                                 , quantile(n_patients_optimal_NO , 0.025) 
                                 , quantile(n_patients_optimal , 0.025))
                    , CI_high = c(quantile(n_patients_notoptimal , 0.975) 
                                 , quantile(n_patients_optimal_NO , 0.975) 
                                 , quantile(n_patients_optimal , 0.975)))
kable(propd)

Method	Median	CI_low	CI_high
Indipendent	15294	14306.450	16671.9
Priority Non Optimal	9402	7803.425	10717.6
Priority Optimal	4476	3918.525	5322.0

Plot as barplot

gpropd <- ggplot(propd, aes(x=Method, y=Median)) + 
    geom_bar(position=position_dodge(), stat="identity" 
             , fill = c("firebrick3" , "deepskyblue4" , "seagreen4") 
            , col=rep("black",3)) +
    geom_errorbar(aes(ymin=CI_low, ymax=CI_high)
                  , width=.2
                  , position=position_dodge(.9)
                  , size=0.2
                  )
gpropd

ggsave(filename="../Figures/fig5C.svg", plot=gpropd, device = "svg")

## Saving 10 x 8 in image

Proportion trial 1-sample

Perform a 1-side 1-sample sample size calculation for difference in proportion with:

• pCase = 0.4
• pControl = 0.1
• power = 0.8

NOT OPTIMAL

set.seed(111)
# NOT optimal
prop3 <- mclapply(1:100 , function(x){
  out <- propPowerSampleSize(PTD_panel
   , alterationType = c("mutations", "copynumber")
   , pCase = 0.4
   , pControl = 0.1
   , power = 0.8
   , type = "arcsin"
   , var = "drug"
   , side = 1
   , tumor.weights = c(brca=brca_tot
                               , luad=luad_tot
                               , lusc=lusc_tot
                               , stad = stad_tot
                               , coadread = coadread_tot
                               )
   , noPlot=TRUE)
  out <- out[ order(out$Var) , ]
  return(out)
} , mc.cores=detectCores())

Screenings <- sapply(prop3 , function(x) x$ScreeningSampleSize)
medianScreening <- apply(Screenings , 1 , function(x) {
  paste(median(x) 
        , paste0( "(" , quantile(x , 0.025) 
                        , "-" 
                        , quantile(x , 0.975) , ")"))
})

prop3DF <- data.frame(Var = prop3[[1]]$Var
                  , ScreeningSampleSize = medianScreening
                  , EligibleSampleSize = prop3[[1]]$EligibleSampleSize
                  , Beta = 0.2
                  , Power = 0.8
                  , Proportion.In.Case.Control = "0.4 - 0.1")
#export
#write.table(surv2, file = "~/Desktop/notoptimal_ptd-showcase.txt", quote=FALSE, sep=",")

#preview
kable(prop3DF[ match(c("Full Design" , drugs , "no_drug") , prop3DF$Var) , ] , row.names=FALSE)

Var	ScreeningSampleSize	EligibleSampleSize	Beta	Power	Proportion.In.Case.Control
Full Design	20 (20-21)	13	0.2	0.8	0.4 - 0.1
HER2 inhibitor	749 (631.275-910.975)	13	0.2	0.8	0.4 - 0.1
PARP inhibitor	147 (135.475-156.525)	13	0.2	0.8	0.4 - 0.1
AKT inhibitor	44 (42.475-46)	13	0.2	0.8	0.4 - 0.1
EGFR inhibitor	526 (458.95-607.1)	13	0.2	0.8	0.4 - 0.1
FGFR inhibitor	132 (125.475-142)	13	0.2	0.8	0.4 - 0.1
ALK inhibitor	382 (343.475-435)	13	0.2	0.8	0.4 - 0.1
BRAF inhibitor	309.5 (281.425-334.525)	13	0.2	0.8	0.4 - 0.1
MET inhibitor	803 (684.225-955.25)	13	0.2	0.8	0.4 - 0.1
NOTCH inhibitor	233 (213.425-255.575)	13	0.2	0.8	0.4 - 0.1
Immune Checkpoint Inhibitor	981 (841-1201.925)	13	0.2	0.8	0.4 - 0.1
no_drug	61 (59-63.525)	13	0.2	0.8	0.4 - 0.1

• Median Number of patients necessary with CI: 4342 (4027.95-4638.175)

OPTIMAL

set.seed(112)
simulprop3 <- mclapply(1:100 , function(x) {
    propPowerSampleSize(PTD_panel
   , alterationType = c("mutations", "copynumber")
   , pCase = 0.4
   , pControl = 0.1
   , power = 0.8
   , var = "drug"
   , type = "arcsin"
   , tumor.weights = c(brca=brca_tot
                               , luad=luad_tot
                               , lusc=lusc_tot
                               , stad = stad_tot
                               , coadread = coadread_tot
                               ) 
        , noPlot=TRUE
        , priority.trial=drugs
        , priority.trial.order="optimal")$`pCase:0.4 | pControl:0.1 | Power:0.8`$Summary
} , mc.cores=detectCores())
prop3summaryMat2 <- simulprop3[[1]]
for(i in rownames(prop3summaryMat2)){
    for(j in colnames(prop3summaryMat2)){
        vec <- sapply(simulprop3 , function(x) x[i , j])
        prop3summaryMat2[i , j] <- paste(median(vec) , 
                paste0( "(" , quantile(vec , 0.025) 
                        , "-" 
                        , quantile(vec , 0.975) , ")"))
    }
}

#save
#write.table(summaryMat2, file = "~/Desktop/optimal_ptd-showcase.txt", quote=FALSE, sep=",")
#preview
kable(prop3summaryMat2)

	Immune Checkpoint Inhibitor	MET inhibitor	HER2 inhibitor	EGFR inhibitor	ALK inhibitor	BRAF inhibitor	NOTCH inhibitor	PARP inhibitor	FGFR inhibitor	AKT inhibitor	Total
Screened	1208 (78.475-1489)	0 (0-1092.775)	0 (0-105.025)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	1233.5 (1053.425-1489)
Eligible	16 (16-18)	18 (16-25)	20 (16-25.525)	29 (23-38)	33 (25.475-40.525)	43 (35-54)	58 (47.475-70)	88.5 (72.425-111.1)	84 (71-104.05)	260.5 (219.9-310.575)	643 (552-777.25)
Not.Eligible	577 (38-707.05)	0 (0-523.675)	0 (0-50.7249999999999)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	0 (0-0)	589.5 (507.325-707.05)

Number of patients necessary under optimal design (median value): 1233.5 (1053.425-1489)

Basically the number of samples screened for Immune check inhibitor, that is targetable only in a very small number of patients, is sufficient to cover all the other drugs in most of the cases.

The percentage of samples saved from passing from a non optimal (all drugs separated) to an optimal design (all drugs evaluated together) is 71.59%

PRIORITY NON OPTIMAL

Perform a priority trial starting from the drug with the highest frequency (NON OPTIMAL)

set.seed(112)
simulprop3_NO <- mclapply(1:100 , function(x) {
    propPowerSampleSize(PTD_panel
   , alterationType = c("mutations", "copynumber")
   , pCase = 0.4
   , pControl = 0.1
   , power = 0.8
   , var = "drug"
   , type = "arcsin"
   , tumor.weights = c(brca=brca_tot
                               , luad=luad_tot
                               , lusc=lusc_tot
                               , stad = stad_tot
                               , coadread = coadread_tot
                               ) 
        , noPlot=TRUE
        , priority.trial=names(drugFreqs)
        , priority.trial.order="as.is")$`pCase:0.4 | pControl:0.1 | Power:0.8`$Summary
} , mc.cores=detectCores())
prop3summaryMat2_NO <- simulprop3_NO[[1]]
for(i in rownames(prop3summaryMat2_NO)){
    for(j in colnames(prop3summaryMat2_NO)){
        vec <- sapply(simulprop3_NO , function(x) x[i , j])
        prop3summaryMat2_NO[i , j] <- paste(median(vec) , 
                paste0( "(" , quantile(vec , 0.025) 
                        , "-" 
                        , quantile(vec , 0.975) , ")"))
    }
}
#save
#write.table(summaryMat3, file = "~/Desktop/optimal_ptd-showcase.txt", quote=FALSE, sep=",")
#preview
kable(prop3summaryMat2_NO)

	AKT inhibitor	FGFR inhibitor	PARP inhibitor	NOTCH inhibitor	BRAF inhibitor	ALK inhibitor	EGFR inhibitor	HER2 inhibitor	MET inhibitor	Immune Checkpoint Inhibitor	Total
Screened	54 (53-56)	124 (114.475-134.525)	66.5 (50.95-83.05)	177 (154.425-218.2)	105.5 (60.95-159.2)	253 (190.175-324.575)	237.5 (118.475-368.25)	610 (425.4-831.4)	257.5 (47-478.05)	638.5 (351.5-965.525)	2534.5 (2231.475-2929.175)
Eligible	723 (635.475-850.075)	177 (153.325-202)	141 (120.475-168.575)	67 (58-82.575)	69 (57.475-83.525)	40 (34-48)	36 (29.475-44)	21 (19-24)	21 (18.475-24.525)	16 (16-16)	1308.5 (1163.8-1521.7)
Not.Eligible	27 (25-28)	60 (54.475-65.525)	32 (24.475-40.525)	85 (74-104.625)	51 (29.475-76.1)	122.5 (92.85-157.05)	115.5 (57.475-181.725)	294 (206.375-396.625)	124 (23-225.925)	630.5 (346.025-954.525)	1532 (1285.85-1847.225)

Number of patients necessary under non optimal priority design (median value with CI): 2534.5 (2231.475-2929.175)

FULL SUMMARY OF THE THREE METHODS

Create a dataframe with median and CI for the three methods

propd1 <- data.frame( Method = c("Indipendent" , "Priority Non Optimal" , "Priority Optimal")
                    , Median = c(median(n_patients_notoptimal) 
                                 , median(n_patients_optimal_NO) 
                                 , median(n_patients_optimal))
                    , CI_low = c(quantile(n_patients_notoptimal , 0.025) 
                                 , quantile(n_patients_optimal_NO , 0.025) 
                                 , quantile(n_patients_optimal , 0.025))
                    , CI_high = c(quantile(n_patients_notoptimal , 0.975) 
                                 , quantile(n_patients_optimal_NO , 0.975) 
                                 , quantile(n_patients_optimal , 0.975)))
kable(propd1)

Method	Median	CI_low	CI_high
Indipendent	4342.0	4027.950	4638.175
Priority Non Optimal	2534.5	2231.475	2929.175
Priority Optimal	1233.5	1053.425	1489.000

Plot as barplot

gpropd1 <- ggplot(propd1, aes(x=Method, y=Median)) + 
    geom_bar(position=position_dodge(), stat="identity" 
             , fill = c("firebrick3" , "deepskyblue4" , "seagreen4") 
            , col=rep("black",3)) +
    geom_errorbar(aes(ymin=CI_low, ymax=CI_high)
                  , width=.2
                  , position=position_dodge(.9)
                  , size=0.2
                  )
gpropd1

ggsave(filename="../Figures/fig5D.svg", plot=gpropd1, device = "svg")

## Saving 10 x 8 in image

One figure for all

sampSizeMethods <- rbind(pfs2 , pfs1 , propd , propd1)
sampSizeMethods$Design <- rep(c("Survival 2 Sample" , "Survival 1 Sample" , "Proportion 2 Sample" , "Proportion 1 Sample") 
                              , c(3,3,3,3))

gsampSizeMethods <- ggplot(sampSizeMethods, aes(x=Design, y=Median , fill=Method)) + 
    geom_bar(position=position_dodge(), stat="identity" , color="black")+
    scale_fill_manual(values=c("firebrick3" , "deepskyblue4" , "seagreen4"))+
    geom_errorbar(aes(ymin=CI_low, ymax=CI_high)
                  , width=.2
                  , position=position_dodge(.9)
                  , size=0.2
                  )
gsampSizeMethods

ggsave(filename="../Figures/fig3B.svg", plot=gsampSizeMethods, device = "svg")

## Saving 10 x 8 in image

Session info

For reproducibility purpose:

sessionInfo()

## R version 3.3.3 (2017-03-06)
## Platform: x86_64-apple-darwin13.4.0 (64-bit)
## Running under: macOS Sierra 10.12.5
## 
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
## 
## attached base packages:
## [1] parallel  stats     graphics  grDevices utils     datasets  methods  
## [8] base     
## 
## other attached packages:
## [1] gdtools_0.1.4                 svglite_1.2.0                
## [3] PrecisionTrialDesigner_0.99.0 ggplot2_2.2.1.9000           
## [5] knitr_1.15.1                  DT_0.2                       
## 
## loaded via a namespace (and not attached):
##  [1] Biobase_2.34.0                httr_1.2.1                   
##  [3] tidyr_0.6.3                   jsonlite_1.5                 
##  [5] AnnotationHub_2.6.5           gtools_3.5.0                 
##  [7] shiny_1.0.3                   assertthat_0.1               
##  [9] LowMACAAnnotation_0.99.3      interactiveDisplayBase_1.12.0
## [11] highr_0.6                     stats4_3.3.3                 
## [13] Rsamtools_1.24.0              yaml_2.1.14                  
## [15] ggrepel_0.6.5                 RSQLite_1.1-2                
## [17] backports_1.0.5               digest_0.6.12                
## [19] GenomicRanges_1.24.3          RColorBrewer_1.1-2           
## [21] XVector_0.12.1                cgdsr_1.2.6                  
## [23] colorspace_1.3-2              htmltools_0.3.6              
## [25] httpuv_1.3.3                  R.oo_1.21.0                  
## [27] plyr_1.8.4                    XML_3.98-1.9                 
## [29] biomaRt_2.30.0                zlibbioc_1.18.0              
## [31] purrr_0.2.2.2                 xtable_1.8-2                 
## [33] scales_0.4.1.9002             gdata_2.17.0                 
## [35] brglm_0.6.1                   BiocParallel_1.6.6           
## [37] tibble_1.3.3                  IRanges_2.6.1                
## [39] SummarizedExperiment_1.2.3    BiocGenerics_0.20.0          
## [41] lazyeval_0.2.0                magrittr_1.5                 
## [43] crayon_1.3.2                  mime_0.5                     
## [45] memoise_1.0.0                 evaluate_0.10                
## [47] R.methodsS3_1.7.1             BiocInstaller_1.24.0         
## [49] tools_3.3.3                   data.table_1.10.4            
## [51] hms_0.3                       stringr_1.2.0                
## [53] googleVis_0.6.2               S4Vectors_0.10.3             
## [55] munsell_0.4.3                 AnnotationDbi_1.36.2         
## [57] Biostrings_2.40.2             GenomeInfoDb_1.8.7           
## [59] profileModel_0.5-9            rlang_0.1.1                  
## [61] grid_3.3.3                    RCurl_1.95-4.8               
## [63] htmlwidgets_0.8               bitops_1.0-6                 
## [65] shinyBS_0.61                  labeling_0.3                 
## [67] rmarkdown_1.4                 testthat_1.0.2               
## [69] gtable_0.2.0                  codetools_0.2-15             
## [71] DBI_0.6                       reshape2_1.4.2               
## [73] R6_2.2.2                      GenomicAlignments_1.8.4      
## [75] dplyr_0.5.0                   rtracklayer_1.32.2           
## [77] rprojroot_1.2                 readr_1.1.0                  
## [79] stringi_1.1.5                 Rcpp_0.12.12

---

1. Frequency is calculated on the the selected population for the study (i.e. frequency sum = 1) and not on the overall cancer population that includes all tumor types.↩

2. Frequency is calculated on the the selected population for the study (i.e. frequency sum = 1) and not on the overall cancer population that includes all tumor types.↩